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Abstract 

The onset of collective behavior in a population of globally coupled oscil¬ 
lators with randomly distributed frequencies is studied for phase dynamical 
models with arbitrary coupling; the effect of a stochastic temporal varia¬ 
tion in the frequencies is also included. The Fokker-Planck equation for the 
coupled Langevin system is reduced to a kinetic equation for the oscillator 
distribution function. Instabilities of the phase-incoherent state are studied 
by center manifold reduction to the amplitude dynamics of the unstable 
modes. Depending on the coupling, the coefficients in the normal form can 
be singular in the limit of weak instability when the diffusive effect of the 
noise is neglected. A detailed analysis of these singularities to all orders 
in the normal form expansion is presented. Physically, the singularities are 
interpreted as predicting an altered scaling of the entrained component near 
the onset of synchronization. 
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1 Introduction 


A phase dynamics description of N limit cycle oscillators is feasible for sufficiently weak 
interaction. ii I a The uncoupled system is attracted to an N-torus and for weak coupling 
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this attracting torus will persist; the oscillator phases ... , 0 jv(O) are coordinates on 

the perturbed torus. In this paper we study such systems where the interaction is mediated 
by a global coupling. BB More precisely, we consider models of the form, 

K N 

di = Wi + -T7 ^2 f (Qj — + &(0 (1) 

7V 3=1 

with the frequencies uof the uncoupled oscillators drawn from a distribution g(ui) char¬ 
acteristic of the population or array. For both physical and mathematical reasons, it is 
interesting to include in (HD the effect of extrinsic white noise £j(t) perturbing the deter¬ 
ministic phase dynamics; this term is defined by the ensemble averages < £i(t) >= 0 and 
< £j(s) £j(t) >= 2 D5ij5(s — t). When the oscillators are uncoupled (.K = 0), then the noise 
and the intrinsic variation in the frequencies spreads the population in phase. For suitably 
chosen coupling functions /(</>) as the strength of the coupling K increases, the instantaneous 
frequencies 6 t may become entrained and the phase disorder of the uncoupled system gives 
way to an entrained state with partial or total phase coherence. 

The collective oscillations produced by such frequency-entrained states is of interest in 
physics, chemistry, biology and rather recently neuroscience. [0’ § 0- 0- 0 The phase dynamics 
model (HD was popularized especially by Kuramoto and his collaborators for the coupling 
/(0) = sin (f) since this choice allowed a mathematically tractable description of the onset 
of entrainment.||] For a coupling /(</>) of arbitrary form, the model is a special case of a 
general normal form obtained by averaging. [[| Very recently, the averaged equations for a 
series array of weakly interacting Josephson junctions has been shown to have the form ([[]) 
with /(0) = sin (f) at leading order in the coupling strength K . || This result extended similar 
studies that had treated the junctions as identical oscillators. || 

Most research on phase oscillators with global couplings has been motivated by their 
appeal as theoretical models of synchronization rather than as realistic descriptions of specific 
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experimental systems. In this literature, until very recently, the preoccupation with the 
simple coupling function considered by Kuramoto has been nearly univeral; see [□]] for 
discussions with many references. Nevertheless, couplings /(</>) of more general form are 
of interest for several reasons. First, interactions that are derived when a reduction to 
phase dynamics is actually carried out can easily have a more complicated structure. 
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Secondly, recent results indicate that when the form of /(0) is generalized, the emergence 
of the entrained state can have qualitative different features. [p~3j Ibj , [l5| More precisely, if 
K c is the critical coupling, above which the synchronized state appears, then the “size” of 
the entrained component scales like (K — K c )9 above threshold. In the absence of noise, 
Daido discovered (5 = 1 when /(</>) contains a second harmonic such as sin 20, but for 
the Kuramoto model the exponent drops to (3 = 1/2, i.e. the entrainment for the more 
general coupling can be much weaker at the same coupling stength /l.|13| In this paper, 
we present a detailed analysis of the first transistion to a state with frequency entrainment 
that illuminates the interplay between the intrinsic disorder g(co), the external noise, and 
the coupling in determining the exponent (3. In particular, we make no further restrictions 
on /(</>) other than to assume it has a Fourier expansion 


OO 

m= £ fne in t (2) 

n=— oo 

For large oscillator populations, it is natural to analyze the system in the limit N — > oo 
because of the relative simplicity afforded by a continuum description of the population using 
a density function p(9,co,t). A convenient definition divides out the frequency dependence 
of the uncoupled system so that g(u>) p(9,uj,t ) d9dco describes the fraction of oscillators with 
natural frequency in [u,u + duj} and phase in [#, 9 + d9 ], and p{9 , a;, t ) d9 describes the fraction 
of oscillators at natural frequency u with phase in [ 9 , 9 + d9]. This latter interpretation fixes 
the normalization / p(9,uj,t) d9 = 1. The dynamics of the system (!|) is formulated as a 
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kinetic equation for p(6,u,t), 


dp , d(pv) _ dV 

dt + 89 89 2 ’ U 

with the deterministic part of the phase velocity expressed as an integral over the popu¬ 
lation 


r2n r oo 

v(6,u,t)=u + K d9' duj'f(9'-9)p{9',uj',t)g{uj'). (4) 

This equation has been applied recently to the Kuramoto model; |T6[-|[T^| our application 
to (|l|) for arbitrary coupling is a generalization of the results for /(0) = sin (j). [17, For 
completeness, the derivation of (^]) - (f|) is outlined in Appendix [E|. 

The kinetic equation has a special solution, 


= (5) 

which describes a population of oscillators spread out unpreferentially in phase. This uniform 
or incoherent distribution is an equilibrium since v(9, u) = uJ+Kf 0 for p = p 0 implies d t p = 0. 
The occurrence of instabilities of p 0 leading to states with nontrivial phase dependence 
is established by a linear stability analysis of p 0 in section [2j and a nonlinear analysis of 
the instabilities in sections 0-0 The nonlinear analysis focuses on the evolution of the 
linearly unstable mode described by the flow on the unstable manifold of the equilibrium; 
the amplitude equation describing this flow are formulated and studied in general terms in 
section [3]. The detailed calculation of the unstable manifold flow given in section [| relies 
on expansions in the amplitude of the unstable mode. The exponent f3 emerges as the 
scaling exponent required to balance the linear instability against the nonlinear terms in the 
amplitude equation. In the regime of weak linear instability, different exponents can arise 
because the behavior of the nonlinear terms varies considerably depending on the form of 
the coupling and strength of the external noise. These features of the expansion are studied 
in sections 0-0 
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Our conclusions regarding the properties of the amplitude equation are valid to all orders 
in the expansion. This level of generality is possible in part because the recursion relations 
for the expansion are considerably simplified by the symmetry of the problem. The model (^) 

- (|j) has SO(2) or 0(2) symmetry depending on details of the population and the coupling. 
The group 0(2) is generated by rotations (3 ■ ( 9,co ) — (9 + (3, u>) and reflections n ■ (0, cu) = 

— (9,lv) which act on functions p(9,u) in the usual way: for any transformation x £ 0(2), 
(x ■ p)(9, cj) = p(y _1 • (9, to)). The dynamics (||) - (|) is equivariant with respect to rotations 
for arbitrary choices of g(u) and /(0); if, in addition, g(u) = g(—u) and /(</>) = 

then the system commutes with the reflection k as well. In the latter circumstance the model 
has 0(2) symmetry, otherwise the rotational symmetry alone corresponds to SO(2). Note 
that the mean frequency of g(co) can always be shifted to zero and we assume that this has 
been done. 

The onset of linear instability for the incoherent state po has some unusual features, 
especially in the limit of weak noise D —> 0. This was first noticed in the context of the 
Kuramoto model, [T7]-[T^] and is emphasized in the linear analysis of section^. In the absence 
of noise, the unstable modes correspond to eigenvalues emerging from a neutral continuum at 
onset; thus, the standard framework of center manifold reduction cannot be directly applied. 
We proceed by first deriving the amplitude equation on the center-unstable manifold with 
D > 0, in this case the diffusion term (S) shifts the continuous spectrum off the imaginary 


axis and the center manifold reduction is straightforward. [pQ , |2lH The resulting amplitude 
equation is then examined by first reducing the noise D —> 0 and then letting the linear 
growth rate of the mode approach zero 7 —> 0 + . In this way, we start with an amplitude 
equation on the unstable manifold and examine its features as the unstable eigenvalues are 
allowed to approach a neutral continuous spectrum on the imaginary axis. 

Our procedure for extracting (3 from the amplitude equation can be illustrated by a 
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simple example. Let A(t) denote the mode amplitude satisfying 


A(t) = 7 A + ai( 7 )A 3 + a 2 ( 7 )A 5 H- ( 6 ) 

where the nonlinear coefficients Gq, a 2 , ... depend on the parameters of the problem such as 
the growth rate 7 . For small 7 , this equation should be viewed as a singular perturbation 
problem with the linear term 7 A representing the perturbation. If 7 > 0, then there is 
always a neighborhood of the equilibrium A = 0 where the perturbation dominates the 
unperturbed system and can completely change the dynamics. As in other such singular 
problems, a possible strategy is to seek a (singular) change of variables which transforms ([|) 
into a regular perturbation problem. Thus it is natural to define a new amplitude by 

A(t) = (7) 

and rewrite (||D in terms of r(r) 

C -!A = 7 1 - 5 r(r) + 7 2 / 3 ~ 5 ai( 7 )r 3 + 7 4 / ^ 5 a 2 ( 7 )r 5 H-. ( 8 ) 

ar 

If possible, the choice of (3 and 5 should be made so that each term in (j 8 j) is well behaved as 
7 —> 0 + , and furthermore so that the effect of 7 > 0 is a regular perturbation of the system 
at 7 = 0. The standard values are of course (3 = 1/2 and 5 = 1 yielding 

dv 

— = r(r) + ai( 7 )r 3 + 7 a 2 ( 7 )r 5 H-. (9) 

dr 

The unperturbed system r = r + Gq(0)r 3 now includes the linear term and near r = 0 the 
effect of small 7 is unimportant. The change of variables (^) is singular at 7 = 0 as expected, 
but nevertheless we have determined the overall scaling (3 = 1/2 for nonlinear solutions near 
the equilibrium in the regime of weak instability. In addition, since terms in (||) at fifth order 
or higher are at least 0( 7 ), the unperturbed problem truncates to a simple balance between 
the linear instability and the dominant nonlinear term. Note that since 7 = 0 at K = K c , 
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then just above threshold 7 ~ (K — K c ), and the definition in ([?) implies |A| rv./ (k - K c y 
for the scaling of the mode ampltude near onset. 

The preceeding discussion contains the unstated assumption that the coefficients op, a 2 , . . . 
have well-defined finite limits as 7 —> 0. In the phase dynamical model (|3;) - (||) this as¬ 
sumption is correct for any choice of /(</>) and g(ui) as long as D > 0. However if we set 
D = 0 before taking the limit 7 —> 0 + , then the nonlinear coefficients are singular in general 
although the details depend on /(0). These singularities require a different choice for (3 to 
accomplish the desiderata enumerated above. It turns out that when one has a coupling 
with only a single Fourier component such as in the Kuramoto model then the expansion 
coefficients are nonsingular to all orders. Thus single component couplings are special and 
very atypical. The most typical or generic behavior is /3 = 1 as was found by Daido in the 

1] We further discuss the relation between our approach and 


examples he considered. [13 
his theory in section [| 

A key aspect of the singular behavior encountered in this problem is the emergence of 
the critical eigenvalues from a neutrally stable continuous spectrum. It has become clear 
recently that this characteristic is shared by a diverse set of additional examples which include 
the onset of linear instability in collisionless plasma |[22l|-|[24|| and in fluid shear flows p5 
26H , instabilities of solitary waves 


27[]-[|29|, and bifurcations in “mean field” descriptions of 


the dynamics of bubble clouds in fluids[|30|. Aside from the present work, only for the 
collisionless plasma models have the resulting singularities in the amplitude expansion been 


fully analyzed. Our proofs in this paper are adapted from the techniques used in [23 
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2 Linear analysis 


Let r/ = p(9,u,t) — p 0 denote the deviation from p 0 , and rewrite the dynamics (|3|) - (|j) to 
describe the evolution of rp. 

= £r) + (10) 

where the linear and nonlinear operators are defined by 

Cr, = C0-(w + A'/„)^ + tt| o 2 '«jO'|” du'fXff (11) 


r2n ro o 

A f{rj) = K / d6' / du'rj(9',uj',t) g(u f ) 

J 0 J —oo 


V(0,w,t) f'{9' - 9) 

p0,u,,t)f(6'-6) 


( 12 ) 


In © - (Q), /'(0) = df /d(f>, and note that the normalization of p implies 


f27T 


d9rj(9, to, t) = 0. 


(13) 


The equilibrium po is invariant under translations [3 and reflections /c; hence the dynamics 
for rj has the same symmetries as the original system. This can be verified explicitly by noting 
that the operators C and A f commute with rotations for arbitrary choices of g(u) and /(0), 
and also commute with n when g{pj) = g(—uj) and /(0) = —/(—</>)• 

Since po is independent of 9 , the eigenfunctions £T = AT have the form 


${9,u) = ip{u)e tne , 


(14) 


with ip(uj) required to satisfy the eigenvalue equation 


Lnlp = Xlp, 


(15) 


where 


L n ip = —in 


/ OO 

du' g(u') ip(u r ) 

-OO 


(16) 







Because the n — 0 component of rj is zero (|T3|), we can restrict our consideration of (JT5|) to 
n 0 . 

The spectrum of L n corresponds to values of A where the resolvent R n ( A) = (A — L n )~ l 
fails to be a bounded operator. By solving the equation (A — L n )(j) = <f for ) we obtain 
the resolvent explicitly, i.e. f = R n (X)f> with 


R n { A)</> = 


—i/n 


</>M - 


Kf: 


du' g(u/) 0 (u/) 


lu + Kfo — iX/n — inD 
where, in the second term, A n (z) denotes 

An 7) = 1 + K f* J_ 


\ n (iX/n ) 7-oo a )' — KJ\) — iX/n — inD 


(17) 


du- 




(18) 


/-oo uj + Kfo — z — inD 
The details of a completely rigorous discussion of R n { A) will depend of course on the choice 
of a specific function space for the problem, e.g. square integrable functions of uj. For our 
purposes, it suffices to ignore these technicalities and discuss the spectrum of L n carefully 
though heuristically. 

There are two ways i? n (A)0 can be singular: when iX/n hits a zero of the function A n 
and when the denominator in (0) vanishes, i.e. when A = —in(uj + Kfo — inD). Since 
—oo < a i < oo, this latter set of values coincides with the line Re(A) = —n 2 D, while the first 
set of values are typically discrete A = —inzj since A n (zj) = 0 will have discrete solutions 
{z\, Z 2 , • • •}• These two components generally correspond to the essential spectrum and the 
point spectrum (eigenvalues), respectively. 

This distinction is clarified by solving the eigenvalue equation directly. Let A = —inz, 


then fll5i) reads 


/ OO 

du’ g(u') if (u 1 ). 

-OO 


(19) 


The possible solutions for Im(z) 7 ^ — nD and Im(£) = — nD are quite different and we discuss 
them separately. 
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For Irn(^) ^ —nD or Re(A) 7 ^ — n 2 D , we can divide by (a; + K /o — z — inD) and obtain 


the form of the eigenfunction 




-K f* n j dtv' g{tv')ift{tv') 
tv + K /o — z — inD 


( 20 ) 


in terms of the constant / dv' gift. A nontrivial solution in ( 20 ) is found only if this constant 
is non-zero, in which case we can adopt the normalization 



du' g{v') ift(tv') 


1 . 


( 21 ) 


Together (|]]) and (21) determine an eigenfunction 

-A' f* 

J n 


lft(u) = 


tv + K / 0 — z — inD 


( 22 ) 


provided the normalization condition (|2l| ) for ift can be satisfied self-consistently. This re¬ 
quirement is precisely A n (z) = A n (iX/n) = 0. Thus a root A n (z) = 0 implies the existence 
of an eigenvalue A = — inz with eigenfunction ([ 2^). 

For all choices of gftv) and f((ft), A n (z) satisfies the identity 


A n (z)* = A-n{z*)] 


(23) 


in addition, with reflection symmetry f((ft) = — /(— (ft) and g{v) = g(—iv), a further identity 
holds 

A n (z) — A_ n (™z). (24) 

The first identity implies that complex eigenvalues must come in conjugate pairs (A, A*) with 
conjugate eigenvectors ^(9,u) and T( 6 I , a;)*, respectively. For a real eigenvalue (Re(^) = 0), 
the first identity reduces to the second (|24[) and implies that the eigenvalue will always 
have (at least) two linearly independent eigenvectors, T and T*. When we have a complex 
eigenvalue for a reflection symmetric population, then (^4j) holds and again there will be two 
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linearly independent eigenvectors ^(0, uo) and T(— 9, —a;). Without reflection symmetry, the 
typical case is a complex eigenvalue, but with reflection symmetry both real and complex 
eigenvalues can occur generically. 


The actual occurrence of roots for A n (z) is a detailed question depending on the frequency 
distribution g and the coupling; examples have been found for which A n (z) has no roots. 13 
In the complex z plane, A n (z) is an analytic function for lm(z) 7 ^ — nD with a discontinuity 
across Irri(z) = — nD whenever g(Re(z) — Kf 0 ) 7 ^ 0. The discontinuity signals a branch 
cut for A n (z) and the magnitude of the discontinuity is given by the boundary values at 
Re(z) = r: 

9 M 


lim A n (r — inD ± ie) = 1 + K ft 

e—>0+ 


PV / dU 1 jy- f 

J —00 00 T A jo 


± in Kf* g{r — Kf 0 ). (25) 


As the parameters of the coupling or frequency distribution g are varied, the roots of A n (z) 
generally change. In particular, roots can appear or disappear at the branch cut. 

When D = 0, the appearance of a root coincides with the onset of linear instability and 
the boundary values ( |25| ) provide closed form expressions for the coupling K c and frequency 
h2 c = — Im(A) at criticality. Let z 0 = (h2 + i^)/n denote the root with 7 > 0, and 0 = 
A n ((h2 + i~j)/n). As 7 —> 0 + , we obtain 


0 


1 + K c f* 



duo 


uo — Q' 


+ inK c f‘g(Sl'). 


(26) 


where = Vt c /n — K c fo- The real and imaginary parts of this equation must vanish 
separately and this implies 


~Min) 

7T l/nIMfT) 


(27) 


ttR e(f n )g(Q,') 

I m (/n) 


jf ^ \g(ti + A) - g(Q' - A)]. 


(28) 
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The second equation determines O' and then the first equation determines /i c ; for n — 1 
and /o = 0 these expressions reduce to the formulas derived by Daido in his order function 
analysis. pi[ 

For Im(z) = —nD or Re(A) = —n 2 D, the eigenvalue equation (0) has no nonsingular 
solutions. One can introduce distributions as solutions following the closely related spectral 
theory for the Vlasov equation, but we shall not require this development for our study (cf. 
|32|| for details). 

In general, the spectrum of £ has both eigenvalues and a continuous component: for 
each mode number n = 1 , 2,..., there is a line of continuous spectrum at Re A = —n 2 D\ in 
addition £ has eigenvalues when the function A n has roots. We are specifically interested 
in eigenvalues that cross the imaginary axis, signifying instabilities of p 0 . Depending on the 
coupling, this can occur for any Fourier component n = l, provided there is a root A; (z 0 ) = 0 
such that Im(zo) becomes positive. Henceforth l and Zq refer to the eigenvalue A = — Uzq 
determining the instability and the unstable eigenfunction is 


V(9,u) = ^(u)e iW 

-Kft 


4>{v) = 


(29) 

(30) 


to + K f§ — zq — ilD 

Clearly the coupling must have a Fourier component fi ^ 0. 

We further restrict the problem with two additional assumptions. First, the root Zq is 
simple, i.e. A[ (z 0 ) 7^ 0 where A[ = dAi/dz] secondly, the center subspace at criticality is 
two-dimensional, in particular, all other roots of A n for any n > 0 have negative imagi¬ 
nary parts and remain bounded away from the real axis. The first assumption ensures that 
there are no generalized eigenvectors; the second restricts the analysis to populations with 
rotation symmetry (SO(2)) and complex eigenvalues, or populations with rotation and re¬ 
flection symmetry (0(2)) and a real eigenvalue. [|Kf[ This restriction excludes one possible 
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codimension-one bifurcation, a Hopf bifurcation with 0(2) symmetry, which involves a four¬ 
dimensional center subspace; a discussion of 0(2) Hopf for the Kuramoto model has been 
given elsewhere. |[TH| 


The projection of rj along the critical mode (BO) requires an inner product 


p2tt roo 

(A,B) = del dcuA(d,cu)*B(9, to), 
J 0 J—o o 


(31) 


so that we can define an adjoint operator for £ 


(£ f A,B) = (A,CB). 


(32) 


In terms of the Fourier expansion A = )T n A n exp(in9), the inner product becomes {A, B) = 
2tt J2 n < A n , B n > using a convenient notation for the a; integration 

/ OO 

duj A n {uj)* B n {uj). (33) 

-OO 

A simple calculation shows that 

OO 

C'A= •£ '"ihM, (34) 

n =—oo 

with 

r /*oo 

(L* A n )(a>) = in (a; + K f 0 + inD)A n (u) + K f n g(u) / du' A n ( u/) . (35) 

L J — O O 

The discussion of the spectrum of L\ parallels the analysis for L n , and we simply state 
the results we shall require. The root A; (z 0 ) = 0 implies an eigenfunction £^T = A*T where 


»( 0 ,W) = —- 0 (c^) e iW 

Z7T 

= 


-9(v) 


(36) 

(37) 


A( (so)*(a; + Kfo ~ 4 + HD )' 

This adjoint eigenfunction satisfies (T, T) = 1 and defines the projection, rj —> (T,?/) T, 
from 77 onto the T component of 77 . This projection is the only aspect of the adjoint theory 
that we need. 
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Before developing the nonlinear theory, we illustrate our discussion with a simple example 
of the instabilities we intend to analyze. Let the oscillator population be described by a 
Lorentzian frequency distribution, 


9{w) = — 


7r Lcu 2 +A 2 J ’ 

then for Im(z) ^ —ID , A; (z) is easily evaluated by contour integration to obtain 

Im {z)<-lD 


(38) 


A i iz) = 


(39) 


For Im(z) < —ID, there is a root at z s = K(f 0 + /*) — i(lD — A) when the coupling satisfies 
K Im(/ { ) > A. The eigenvalues A = — ilz s have negative real parts and fall to the left of the 
continuous spectrum at Re A = —l 2 D. An eigenvalue sitting to the right of the continuous 
spectrum requires a root with Im(z) > — ID ; for this example A; (z) has one such root at 
z u = K(fo + /*) — i(lD + A) when the coupling satisfies Klm(fi) < —A. This root yields 
the eigenvalue 

A = — ilz u = 7 — iQ, (40) 


where 


7 = —l [ID + A + K Im(/;)] (41) 

= IKlReift) + f 0 ]. (42) 


These modes are linearly stable for ID + A + Klm(fi) > 0 and become linearly unstable 
for K > K c where K c = — (ID + A)/Im(/;). Note, that we find eigenvalues of Li only when 
|iFIm(//)| > A; in the limit \K Im(/;)j —> A these eigenvalues approach the continuous 
spectrum and are “absorbed” as the underlying root z u or z s crosses the branch cut onto a 
different Riemann sheet of A^ (z). 
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The diffusion due to external noise simplifies the spectrum by displacing the lines of 
continuous spectrum into the left half plane. With D > 0 the eigenvalues associated with 
instabilities first emerge from the continuum and then subsequently reach the imaginary 
axis. Thus, the critical eigenvalues are isolated on the imaginary axis, and one has a critical 
spectrum of standard form for bifurcation theory. If we let D —> 0, then the unstable 
eigenvalues emerge directly from a continuous spectrum on the imaginary axis. One of our 
central concerns is the effect of this continuous spectrum for D = 0 on the 7 —> 0 + limit of 
the vector field on the center-unstable manifold. This is investigated below by deriving this 
vector field first for D > 0 and then analyzing the effect of taking the limit 7 —> 0 + when 
D = 0. 

3 Amplitude equation on the unstable manifold 


With diffusion included, the appearance of an unstable mode in the spectrum can be analyzed 
by a straightforward application of center manfold reduction. |20|. ^1| This enables us to derive 
approximately the vector field on the center-unstable manifold of p 0l and study the resulting 
bifurcation. This vector field will depend on the parameters of the problem and we shall 
examine its features as the growth rate of the unstable mode is progressively weakened, a 
limit which we denote by 7 —> 0 + . The properties of this limit are strongly dependent on the 
form of the coupling /(0) and whether the external noise D > 0 is held fixed as 7 vanishes. 

Center manifold reduction yields an amplitude equation describing the time-asymptotic 
behavior of the unstable mode. The mode amplitude a{t) = (\k,?/) and describes the pro¬ 
jection of?/ onto the unstable eigenspace (|3C|): 


rj(9, u, t ) = [a(t)T( 6 l , u) + c.c.] + S(0, u;, t) 


(43) 
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where (T, S) = 0. In terms of a and S', the evolution equation (I0[) becomes 


a 

dS_ 

~di 


Xa+ (T, A/"(r/)) 


CS + Af(v) 


T + c.c. 


(44) 

(45) 


These coupled equations are equivalent to (01; by restricting them to solutions on the 
center-unstable manifold we obtain an autonomous two-dimensional flow for a(t). 

Near the equilibrium rj = 0, the center-unstable manifold is described by a function H , 


S(6,uj) = 


(46) 


that satisfies 0 = H(9, to, 0, 0) = d a H(9, u;, 0, 0) = d a * H(9, to, 0, 0). Using this representation 
of the manifold, the general solution in ([43|) can be adapted for solutions on the manifold, 

r] u (6, u, t ) = [a(f)T(6 l , u) + c.c.] + H(9, u, a(t), a(t)*); (47) 


the time dependence of i] u (9,uj,t) is entirely determined by the dynamics of a{t). An au¬ 
tonomous equation for a(t) follows next by restricting the general equation ( pH| ) to solutions 
of the special form (06 

d = Ao+(f,jV(i|“)). (48) 


This is the amplitude equation we wish to derive and analyze. 

Certain general features of the amplitude equation follow from the symmetry properties 
of the evolution equation ( |TQ| ) which are, of course, shared by its expression in (|44]) - (|45|) . 
The action of the rotations (3 and reflection k on the variables (ex, S(9, u)) can be determined 
from their definitions in (|43|). The rotation f3 ■ ( 9,u ) —> (9 + /3,u) of the phase acts by 


p ■ (a, S(0, ui)) = (a e-y S(0 - f3, u)), 


(49) 


and 

k ■ (a, S(9, u)) = (a*, S(—9, —a;)) (50) 
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for the reflection. 


The right hand side of the amplitude equation commutes with these transformations of a 
whenever they are symmetries of i [Tl] i - (fK|); in particular the amplitude equation commutes 
with rotations (fIU|). There are useful standard results on the expression of such equivariant 
vector fields. For example, a two-dimensional vector field a = V(a, a*) that commutes with 
a —> ae~ ll/3 can be written as V(a,a*) = ap(a) where a = |ct| J and p(a) is a smooth 


function. [33] Hence we know the amplitude equation (fig) has the form 


Aa + (T,A f{j] u )) = ap(a) 


(51) 


although p(a) must be determined from the model. In general, p(a) is complex, but when 
the model has reflection symmetry then p(cr) must be real-valued. 

3.1 Analysis of H(6, a;, a, a*) 

The dynamical invariance of the center-unstable manifold implies an equation for the graph 
function H. Consistency between the time dependence of S in (|46j) and the evolution of S 
described by (fh|)-(|45|) requires that H(9,uj,a(t),a*(t )) satisfy 


du mi 

“&r + “ 9^ 


= CH + N'(r) u ) - [(T, M{r] u )) T + (T, A/'(t 7 u ))* T* 


(52) 


where a on the left hand side refers to the amplitude equation (fi8|). 


There are general constraints on H due to the symmetries of ( [44| ) - (|45|) because we 
expect the unstable manifold to be mapped onto itself by a symmetry transformation. For 
a rotation /3, this symmetry invariance of the manifold implies 


H{6 - u, a, a*) = H{8, u, a e~ up , a* e il(3 ), 


and if reflection symmetry holds then H must satisfy 


(53) 


H{—9 , —a;, ck, a*) = H(9, u, a*, a) 


(54) 
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as well. 


These relations constrain the form of the Fourier expansion of H 

OO 

H(6,cu,a,a*) = J2 H n (co,a,a*) e ind ; (55) 

n =—oo 

applying ([53]) to the components H n shows that the component must vanish unless n is an 
integer multiple of l. Note that the n — 0 component of H is also identically zero because 
of the general property / ddg = 0. When n = ml > 0, then H m i must have the form 

Hi(co,a,a*) = aahi(co,a) (56) 

H m i(uj,a,a*) = a m h m (uj,a) m > 1. (57) 

The functions h m (iu,a) are unconstrained by the rotations, but must satisfy h m (—uj,cr)* = 
h m (u),a) when reflection symmetry holds. In addition, the general relation (T, S) = 0 
becomes (T,i7) = 0 for solutions on the unstable manifold, and for the components of H 
this implies 

< ip, h >= 0. (58) 

The symmetry argument for ( |56| ) and (0) can be briefly summarized. Applying (|53| ) to 
the Fourier coefficients of H leads to the identity 

H n {u, a, a*)e~ in(3 = H n (uj, ae~ ilf} , a*e ilf3 ). (59) 


For /3 = 2n/l this immediately implies that H n = 0 unless exp(— i2Tm/l) = 1, hence n must 
be an integer multiple of /. Recall that if a function F(uj , a, a*) is rotation invariant under 
phase shifts a —> a exp then it can be written as F(u), a, a*) = g(u), a) where the form 
of g(co,a) is determined by -F.|33|] In this case, the invariant function (a*) m H m i(u>,a,a*) = 
g{uj, a) must have the form g(co, a) = a m h(uj, a) to accomodate the factor ( a*) m . Dividing by 
( a*) m gives a, a*) = a m h(LU,a ) as in fl57|). The general property d a F[(0, u, 0, 0) = 0 

requires that h(co,a) for Fli have the special form h{u>,a) = ahi(oj,a) as in 
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The equations for h m (us, a) follow by substituting (|56| ) - ( |57p and (|5T| ) into the equation 
for H (^). The results of this straightforward calculation for m — 1,2 and m > 2 are 
summarized below. In writing these equations two further notations prove convenient. First, 
define T m by 

/ OO 

du g(u)h m (u,a), (60) 

-OO 

and secondly let Px denote the projection operator 

(P_l0)M = (t>{u)~ < -0, 0 > i>(u) (61) 

which projects 0 onto the orthogonal complement of the eigenvector 0. 

The n = l component (m = 1) of (^) can be written 

dh f 

[U ~ 2p-p*}h 1 - (p + p*)a-^ = (: 2mKl ) P ± l //(1 + aY\)h 2 + / 2 *)r 2 (0* + crh\) 


+S 21 oT 2 h 2 -\- S','rt , 


(62) 


■5> 

m =3 


m— 1 


fml^mhm.+l + f(m+l)l^m+lK 


The n = 21 component of ( |52|) is 

dh 2 


[L 2l - 2p\h 2 - (p + p*)a 


( 2niK(2l )) | /*[0 + a (hi + Ti0) + n 2 T 1 /i 1 ] 
+fi a ( 1 + ^r*)/i 3 + /3)cr[r 3 0* + oT 3 /i*] 

OO OO 

E ° m fmi r t'w2 + E 


(63) 


* 

m—2 


m=2 


are 


ra=4 


'm+1 


For m > 3, the Fourier components of 

[L mi - mp}h m - (p + p*)cx-^ = (2 tt iK(ml)) < f*[hm-i + oTih m _i] + //cr(l + <jT\)h 

m—2 

+ f(m-n)l^m-nh n ( 64 ) 

n=2 

+/(m-l)*rm-l[0 + O^l] + /(m+l)i a r m+ i[0* + Cr/i*] 

OO OO 'j 

E + E / ( ;, + „)^”r m+ ^; • 


n=2 


n—2 
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Here and below, it is understood that sums are omitted when the lower value of the sum¬ 
mation index exceeds the upper value, e.g. for m — 3 the first sum in ( |64| ) is omitted. 

The symmetry constraints imply that solutions on the unstable manifold ( 0 ) have the 
general form 


p u (0,uj,t) = 


a;) + a(t) a h\(uj, a(t)) e lW 

OO 

+ Y. a m (t)h m (u,a(t))e imW 


(65) 


m =2 


+ c.c. 


with a{t) satisfying the amplitude equation (|48|) . The set of equations (|5T]) and (^2|) - ( |6d|) 
determine the functions p(cr) and h m (u,a), however they cannot be solved explicitly except 
using the amplitude expansions introduced in the next section. What has been achieved to 
this point is to formulate the problem entirely in terms of functions of a single real variable 
cr. This simplification is essential to the analysis to follow. 


4 Expansions and recursion relations 

Our conclusions regarding the nonlinear development of the instability are based on ampli¬ 
tude expansions for p and h m {uj , cr), and the corresponding expansion of T m from (|60|): 


OO 


p(<?) = 

(66) 

3=0 


OO 

h m (uj,a) = u)a 3 

(67) 

3=0 


OO 


r m (a) ^ 9 ■) h"m ^ ^ r m ,j ^ 

(68) 

j =0 



where T mj - = f dujgh m j. The orthogonality condition (|58|) implies 

< ip, hij >= 0 (69) 

at each order in the expansion of h\. 
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The coefficients pj and h m j(u>) are determined by recursively solving (|5T|) and ( |6^ ) - (|64p. 
Relations for pj follow from ( |51|) using the general solution form in (|65| ) and the expansions 
(pb|) - (p%[). At each order in cr\ j > 1, we solve (pT|) for pj and find 


Pj = -(2niKl) < fi 


+/: 


21 


j~ 2 

< 0, h 2J -1 > +J2 r *,n < Vb h 2,j-n-2 > 
n =0 

3-2 

T 2 J-1 < 0,0* > + E F 2,n < 0, K,j-n -2 > 


(70) 


77 -— 0 


3-2 


+ f'U E r^ n < 0, h 3 j - n - 2 > 


77 .— 0 


j+l j-m+1 


j-m 


E /;, E + E fml E r m,n < IP: h. 


m-\-l,j—n—m ^ 


m =3 77=0 777=3 77=0 

The functions h m j(uj) are determined similarly by substituting (|66|) - (p^) into the com¬ 
ponent equations (p2|) - (|64[) . At every order, h m _j is then obtained by applying the resolvent 
operator © to certain auxilliary functions constructed from lower order terms, i.e. 


^ 777 ,j f^ml (fJ'mJ ) Im,j (^) 


(71) 


where 


Pm,j = (nr + j)X + jX* + <5 m ,i(A + A*), 

and the detailed expressions for I m j are provided below. 

The action of the resolvent is given in (p~7|) 


(72) 


A/.// (Pin.] ) Pll.j 


—i /ml 


to + Kf 0 — Vm,j — imlD 


Im,j (^) 


KF ml r du/g(oS)Imj(u/) 


A7-oo u' + Kf 0 - u mJ - imlD 

(73) 


with = i /i m .j /ml. ft is important to note from (|72| ) that z/ mj - can be rewritten 

2 i(j + <5m, 1 ) 7 . 


— A) + 


ml 


(74) 


hence for m > 1, the poles u = — K f 0 + imlD in ( |75| ) all sit in the upper half of the 

complex u plane along the line Re(u;) = Re(z 0 ) — K f 0 . 
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The coefficients p 3 are determined by the integrals T m ,j, < Vb Kn.j >, < bb hm. 3 >, and 
< ip, ip* >. As we are interested in the form of these coefficients near the onset of the 
instability (7 —> 0) when the noise is neglected (i.e. D — 0), it is important to determine the 
behavior of such integrals in this limit. In particular, we will need to estimate any singular 
behavior exhibited by each type of integral. 

In the remainder of the paper, our analysis of the amplitude expansion will be facilitated 
by two identities. From the definition of T m j in (^) and the resolvent expression for h m j(u) 
in @). the following useful relations are obtained 

= *7^ f * 7 + ,^ )W A m (75) 

A ml\^m,j) J ~°° T A Jq Vnn,j tWfllD ) 


hm.j (^) 


'i'lrn. j (^) 


Kfml r m,j 


ml(uj + Kf 0 — — imlD) (uj + Kf 0 — u m j — imlD ) 


(76) 


Given / mj -, these identities allow the properties of h m j and T m j to be quickly determined. 

The I m j are in turn calculated from the following recursion relations in terms of lower 
order coefficients. The application of these recursion relations to calculate the coefficients pj 
is summarized in Table 1 and illustrated in Section || 


4.1 Recursion relations for / 


m j 


The expansion of (|62| ) using (^ 6 ( ) - (| 68 |) determines hij at each order in 
—Ri(pij) Iij where for m = 1 and j > 0 , 


3- 1 


h,j = K 2 + k )Pj-k + (1 + k )p*_ k ] h 1>k 


■ k 1 j ~~ 


(77) 


k =0 


+{2mKl) P ± /, 


j -1 

k 2,j + X! ^n^lj-n-l 

n =0 


+/: 


21 


3 - 1 

^^2 ,j + K,rX2,j-n-l 

n =0 


3 ~ 1 

+ .f '21 k 3,n^l,j-n-l 

n =0 
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jf+2 j— m+2 jf+l j— m+1 

“I" E/ fml E/ ^m—l,nhrn,j—m—n+2 + E/ fml E/ ^m+l,nr m j_ m _ n _|_i 

m=3 n=0 m=3 n=0 

Since < -0, Pj_0 >= 0 by construction, this recursion relation together with the orthogonality 
property of h i ^ (|55|) implies a corresponding orthogonality property for 


< -0, Jij >= 0. 


(78) 


For m = 2, the expansion of (|63D yields h 2 j = —Ri{p 2 ,j) h,j where for j = 0 


h,o = 2mK{2l) f* 


(79) 


and for j > 1, 


j - 1 


= E [ 2 Pj-fc + + Pj-fc)] h 2,k 


k =0 


+2niK(2l) { ft 


+fi 

+hi 


j~2 

h,j-l + E ^3,n rtj_ n _ 2 

ti=0 

3-2 


h,j-i + ^l,n r 1 j— n —2 


71—0 


i-2 


(80) 


+ hi E ^4,n 


71—0 


i-2 


V ; *r 3j _ 1 + 5Z K,n^3,j-n-2 


n =0 


J-3 


+ hi E h,n ^tj-n-2 


71—0 
i 3~k 


j +2 j-fc+2 

+ £/« E h*k— 2,n ^-'k,j—n—k+2 "7 5 ) fkl 5 ) ^fc+2,n 

fc=4 n=0 fc=4 n=0 

For m > 3 and j > 0 we have h m j = 7nj with 

j'-i 

E 

k=0 


Im,j 5 ) \hllpj _/,, “I - k(pj_k hm,k 


r 

J-2 


/* 

hm+lj —1 "F 5 ) ,j—k—2 

+/r 

l 

k =0 



J-l 

m—2 


3~ 1 

hm—l,j "F 5 ) tlm— ljfcbi ,j—k 
fc =0 


+ /(m- 


1)» 


"0r m _4j -f" 5) m-lj'-fc-l 
fc =0 


+ 5Z fnl E/ h'm—n,kT'nJ—k 
n=2 k =0 


(81) 
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+/( 




3 -2 

^T m+ ij_i + ^ /i^ fc r m+ ij_fc_2 
/c =0 

i j-n 


m+j j—n+m 

+ f nl ^m+n^nj-n-k + fnl ^n-m^n,j-n-k+m 

n =2 fc=0 ?i=rrt+2 fc=0 

In this final expression it is understood that a term with a negative subscript is omitted, e.g. 
the term ft /i m +i j_i is dropped when j = 0 . 


5 Pinching singularities at low order 


At this point it is instructive to evaluate the low order coefficients in the expansion and 
understand how the singularities arise mathematically. We first calculate the cubic coefficient 
Pi and show that it has a singularity when f 2 i ^ 0. Then we evaluate the fifth order coefficient 
for the special case of a coupling with f 2 i = 0 and no singularity in pi. 


5.1 Evaluation of the coefficient pi 

The calculation of pi in (Jj7(|) yields 


Pi = -2mKl ft < h 2) o > +/ 2 *,r 2 , 0 < Vh " 0 * > 


(82) 


where the function h 2 0 is determined from ([76]) and (|79|) 

*^2,0 A f*i r 2 , 0 


h 2 ,o{u) = 


2l(u> + Kf 0 — zq — i2lD) (a; + Kf 0 — zq — i2lD) 


(83) 


with r 2j0 given by 


r 2 ,o — 


i/ 2 l [ g{u) 


A 21 (zo) 


(a; + Kf 0 — zq — i2lD) 


2k K ft (z 0 ) 
A 21 {zo) 


+ 0(D). 


(84) 


The asymptotic expression for small D in (|84|) shows that generically T 2)0 7 ^ 0. 

The most important qualitative feature of / 2 > 0 (oj) and h 2 ,o( w ) is that they are meromor- 
pliic functions on the complex u plane. Let Zq = (12 + i^)/l denote the real and imaginary 
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parts of the root, then each function has a pole at u> + = {VL/l — Kf 0 ) + i{l 2 D + 7)//. In 
addition, /i 2 ,o(^) has a second pole at (0/Z — A'/o) + i(2l 2 D + 7 )/l. Both poles lie in the 
upper half plane and although they can approach the real axis as the 7 —> 0 + , the integrals 
< 0, h 2j o > and r 2>0 in (^) have well-behaved nonsingular limits even if D is zero. 

This observation regarding analytic structure can be extended considerably: all of the 
functions and h m j( 00 ) are meromorphic functions in the complex u plane and all of 

their poles he along the vertical line Re (a;) = — Kf 0 + Re(z 0 ). When D — 0, all of the poles 
approach the real axis in the limit 7 —> 0 -1 "; in general, however the poles will be found in both 
the upper and lower half plane and integrals over or h m j(uj) will exhibit pinching 

singularities. These conclusions follow from the recursion relations in light of two simple 
observations. First, each I m j(uj) is built up as a sum of various h m j(u) from lower order 
and the sum of two meromorphic functions is again meromorphic. Secondly, the application 
of the resolvent ([73]) to a meromorphic yields a meromorphic function for h m} j(uj) in 

m- 

The simplest pinching singularity is illustrated by < r ip, r ip* > in pi ; from (|3Dj) and (0) 
one sees that the integrand has poles uj± = (Q/l — Kf 0 ) ± i(l 2 D + 7 )/l above and below the 
contour along the real axis. Now the 7 —> 0 + limit produces a pinching singularity when D 
is zero, and this pinch contributes a factor (7 + l 2 D ) to the denominator when the integral 
is evaluated: 


< - 0 , - 0 * > = 


ifiMfi) 


(85) 


(7 + pb)|/,pa;( Z o )' 

This singularity will occur in the p\ coefficient ( |52| ) provided Dm 0 and / 2 / 7^ 0. 

The derivation of (]S5|) is simple and illustrates the key idea used below in the general 
analysis: make a partial fraction expansion of the integrand to isolate the singularity. In 
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this case we have 


< ip, ip* >= 


Kfi 


du ■ 


9(u) 


A[ (z 0 ) 7-oo (a; — u + )(u — a>_) 
The partial fraction expansion yields 

1 -il/ 2 / 1 1 


( 86 ) 


(87) 


(a; — a; + )(a; — u_) 7 + l 2 D \ (a; — u + ) (u> — u>J)) 1 

so that ( |HB| ) becomes 

(A,(z„)-1 A,w-1\ 

<*’* >_ 2A [(zMy + PD) [-HJr - KIT) ' (88) 

This reduces to (|85|) since A/(^ 0 ) = 0. For D = 0, the cubic coefficient diverges like pi ~ y -1 
as 7 —> 0 + provided / 2 z A 0; this turns out to reflect a general pattern pj ~ y -2 - 74 " 1 which is 
proved below in Theorem | 6 . 1| . 

The significance of this divergence can be appreciated by considering the amplitude equa¬ 
tion truncated at third order 

a = Xa + pi a\a\ 2 , (89) 

and re-examining the determination of (3 in (^). Introducing polar variables a = Apt) exp(— iippt)), 
the amplitude equation reads 


A = 7 A + ReQty)/! 3 


ip = — Im(p!)A 2 


(90) 


where A = 7 — ifl. Let p\ = cqy v denote the asymptotic behavior of the cubic coefficient 
with v > 0 , then in terms of the rescaled amplitude ( 0 ) the equations become 

ClT* 

— = 7 l ~ 5 r + 7 2 ^ -l/-< 5 Re(ai)r 3 ip = hi — y 2 / 3 _ 1 Tm(ai)r 2 . (91) 

ar 

In order to remove the singular perturbation due to the linear term we again take 6 = 1 and 
then solve for (3 — (1 + i /)/2 to balance the cubic nonlinearity against the linear instability. 
This makes it clear that a divergence in the cubic coefficient will increase f3 and for v = 1 we 
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predict (3=1. This of course neglects possible effects from higher order nonlinear terms, in 
particular singularities at higher order could turn out to be more important than the cubic 
singularity. The effects of higher order singularities are systematically studied in section |B|. 


5.2 The fifth order coefficient p 2 


When f 2 i = 0 the cubic coefficient is nonsingular, and we briefly consider what occurs at 
fifth order in the amplitude expansion. Setting f 2 i = 0, the fifth order coefficient is given by 


P 2 = 


(2mKl) {fi [< 0, h 2 , i > +T* h0 < 0, h 2 ,o >] + < 0, h* 2 G >} 


where 


I i,oM = {2iriKl) Pi {/;/i 2 ,o} — (2niKl) fi h 2l o >} 

hAA = 2niK(3l) h 2 AA 

h,iA) = 2pi5 2 ,oM + 2niK(2l) [fih 3i0 (u) + fAAoA) + UoV'M) 

+/jr . 


(92) 


(93) 

(94) 

(95) 


There are five integrals to consider in (|92|) . Since h 20 has poles only in the upper plane, 
the integral < 0, h 2 $ > is not singular. From ([h|), we note that I^.q similarly has poles only 
in the upper half plane hence r 3i0 has no pinching singularity. In addition, / 10 in (^) has 
poles in the upper plane only and this implies that r 1>0 is nonsingular, but the calculation 
is more subtle. From (|75p we have 

r - i / 1 [°° A A,oM 

1,0 A0^i, o ) J-oo (u + Kf 0 - v i i0 - UD ) ' 

Since Aj(z/ li0 ) = A;(z 0 + 2 7//) from Eq. O, the first factor in (^) is singular as 7 —> 0 + . 
However the integral is nonsingular in this limit and because of the orthogonality relation 
(|78|) the leading term is zero so the integral is actually 0( 7 ). The vanishing of the integral 
as 7 — > 0 + cancels the singularity from A;(zq i0 ) and is well-behaved. 
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The remaining integrals do exhibit pinching singularities which can be evaluated with 
the partial fraction procedure illustrated above. Rewriting p 2 to isolate these singular terms 
we find for D = 0: 

p 2 = —(2 iriKl) | fi < x /5,/i 2 ,i > +1^3,0 < Vs ^ 2,0 > + nonsingular terms j (97) 

with 

<tp,h 2} i> = -2 7 rR/*r 3: o / dui(uj) i(u) nonsingular terms (98) 

7-00 u + A /o - (20 + 17/t) 


<vU* 0 > = —2nKfi r 


(99) 


1-00 U + K /o — Zq 

A partial fraction expansion of each of these integrals allows an explicit evaluation of the 
singular behavior 

r°° du y(ujy y(uj)* _ u 2 fi (im(/j) + 0(7)] 


1-00 UJ + A /o — (zo + * 7/0 


3 \ fi\ 2 A' l (z 0 ) 7 2 


( 100 ) 


/■“ dw4>(u})* -0M* = [Im(/i) + 0 ( 7 )] 

7-00 u + Jl /o - So 2 I fi 1 2 A) ( 20 ) 7 2 

each integral diverges like y~ 2 . 

The resulting evaluation of p 2 yields 

21 (tt l K fi) 2 r 3;0 [Im(/Q + 0( 7 )] 

P2 3|/,pa;c- 0 )V 

finally, the T 3j o factor is generically nonzero and can be evaluated from ([75]) as 


+ nonsingular terms; 


r 3 ,o= -2(7r/f/,*) 


*,2 a;'(*o) 




( 101 ) 


( 102 ) 


(103) 


Hence p 2 ~ 7 2 for the special case f 2 i — 0 when p 3 is nonsingular. However, if f 2 i and / 3 ; 
both vanish, then the divergent terms in (| 102 |) disappear and both p\ and p 2 are nonsingular. 











These conclusions are generalized below. In Proposition Q we show that when / 2 ; = 0 
but fsi 7 ^ 0, thenpj ~ y -2 - 7 ' 1 ' 2 describes the strongest possible singularity of pj. In Proposition 
0 we hnd that the coefficients pj are nonsingular at every order if f n i = 0 for all n > 1 . 
The implications of the p 2 divergence for the scaling of the amplitude, 7 ^r( 7 i), can be 


seen by rewriting the radial equation (|91|) to allow for both the third and fifth order terms: 


= 'jr + Re(pi)7 2/3 V 3 + Re(p 2 )7 4/3 1 r 3 . 

dr 


(104) 


When pi is nonsingular and p 2 ~ 7 2 , then we must take (3 = 3/4 to obtain an expansion 
that is finite through fifth order terms. In this case however the general singularity result 


Pj ~ 7 


—2J+2 • 


in Proposition o does not exclude higher order terms that would require a still 


larger value of (3\ this point is discussed further in Section 7.2 below. 


6 Singularity structure of the expansion to all orders 


Our goal is a systematic understanding of the singularities in pj to all orders in the amplitude 
expansion. These singularities are smoothed by the noise (cf. the singularity in (^)), and 
in the discussion below we always set 0 = 0. In this analysis, it is useful to introduce an 
index which keeps track of the worst case singularity for each coefficient. 

6.1 Definition of the index 


The pinching singularities at every order occur in integrands of a standard form. Let 

1 


D n (a, u) = 


(105) 


(a; - aq) (a; - a 2 ) ■ ■ • (a> — a n ) 

for n > 0 where a = (aq,..., a n ) and define D 0 (a,iv) = 1. Evaluating pj for j > 1 involves 
integrals of the form / doj (ft(u)Q(uj) where 


Qipj) = D m (/3,uj)* D n (a,u) 


(106) 


29 





with m + n> 1. Here 0(cn) denotes a smooth function assumed to be well behaved as 7 —> 0 
and to vanish for |u;| —> 00. For example, 0(cn) = g{ uj) is often the smooth function of 
interest. The poles in (|106|) are given by 


Zo - Kfo + i'ySj 

j = !,-• 

• , n 

(107) 

Zq ~ Kfo - i-fCj 

J — 1) • ■ 

. ,m; 

(108) 


hence they lie along the vertical line uj = Re(zo) — Kfo at locations determined by the 
numbers Sj > 0 and Q > 0 which are assumed to be independent of 7 = ilm(z 0 ). 

The index of Q( uj) in (|106|) is defined to be 


Ind [G] = m + n — 1 . 


(109) 


We extend this definition to assign an index if G(uj) is multiplied by a smooth function <f>(uj ) 
as described above 

Ind [<fr(u)g(u)\ = Ind [G{pj)\. (HO) 

Since we assume rri + n > 1 in ( |106| ), Ind [<j>G\ > 0 . If mn 7 ^ 0 , then as 7 —> 0 + , the 
integral / du> <fr{u)G (oj) may diverge due to a pinching singularity at uj — Re(z 0 ) — Kfo■ The 
smooth functions 4>(u) are always assumed to be sufficiently well behaved at large \oj\ that 
the only possible divergence is due to the pinching singularity. The index of <f>G is simply 
related to the maximum possible strength of this divergence. 


Proposition 6.1 ForQ(oj ) in ( |106| ) withm+n > 1 and<f>{ uj) a smooth function as described 
above, the integral of <f>G satisfies 


lim 

7—>o+ 


7 


du (p(u>) Gfjj) 


< 00 


( 111 ) 


with J = Ind [ 4 >G]- If J is replaced by J — 1 ? then the limit (| 111 |) diverges in general unless 
mn = 0 in which case the limit is zero for any J > 0. 
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Proof. 


See Appendix |A]. □ 

Our application of Prop. O to the recursion relations for I m j requires a further gener¬ 
alization of the index to allow for sums of functions with well-defined indices and products 
of (pQ with singular functions of 7 . In each case, the generalized index is defined so that 
Eq. (|1 1 1|) remains true, i.e. the index for the composite function gives the maximal possible 
divergence of its integral. 

First, let Gi(u;) and G 2 (u;) denote functions with well-defined indices, Ind [Gi] > Ind [G 2 ]. 
The index of their sum is defined as the largest of the two individual indices: 

Ind [Gi + G 2 ] =Ind [Gi]. ( 112 ) 

Clearly, Ind [Gi] gives the maximal possible divergence of / du (Gi + G 2 ). Secondly, let g(y) 
denote a function of 7 with the asymptotic behavior 


lim q( 7 ) ~ 7 17 , 

7—>0+ 

then the index of q( 7 ) Gi(u;) is defined by 

Ind [qGi] = Ind [Gi] + v. 


(113) 


(114) 


This completes the definition of the index. 

By applying (|112|) and (|114|) the indices of I rnj and h m j are determined from the recursion 


relations; some simple examples are provided in section |67|. We stress again that when 
applied to a composite function G(co) the estimate in ( | 111 |) does not necessarily determine 
the true singularity, but only an upper bound on the possible divergence of the integral. In 
many examples considered below this upper bound is actually realized, but not always. 
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There are two immediate implications of the index definition worth noting. First, complex 
conjugation doesn’t change the index 

Ind [G] = Ind [G *]. (115) 

Secondly, if G(u>) has a well defined index, then dividing G by (a; — a) or (uj — (3*) simply 
increases the index of G by one: 


Ind [G/(u - a)} = Ind [G/(u - (3*)] = Ind [ G\ + 1. 


(116) 


Here a, (3 introduce additional poles as defined in (|107|) - (|108|) and thereby increase the 


index of every term in G by one; this implies (116). 

The second observation ( |116| ) provides a useful relation between the divergence of T mj - 
and the divergence of < ip, h m j >. Suppose that 7 _iy , v > 0, denotes the actual divergence 
of r m j. Then, from (|75| ) this must be due to certain terms in in addition when m — 1 
then is of order 7 for small 7, which increases the divergence by one. From (El), 

< -0, hmj > is given by 

du> g(u>) 


F ; hrn.'j 


d J- 


ml A[(z 0 ) 7-oo (a; + Kf 0 - z 0 )(uj + Kf 0 - is mJ ) 


(117) 


K fml 


du g(co) 


A'i(z 0 ) 7-oo (a; + Kfo ~ *o)(^ + Kfo ~ u m,j) 

and we conclude that its divergence cannot be worse than y _ (^+ 1 ^ 5m ’ 1 ). The second integral 
in (|1 1 7|) has no pinching singularity and the singularity of the second term is determined by 




7 v . The first term in ( |117|) involves an integral whose integrand differs from the 


integrand in ( |76D only by the extra factor (oo + Kfo — Zo) in the denominator. From ( |1 16|) this 
factor simply increases the index of every term in / (uj + Kf 0 — In particular, 

the terms that determined the singularity of T mj - ~ y -1 " in the first place have their index 
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increased by 1. Thus we can infer the divergence of the first term in (|117| ) by adding one 
to v and subtracting off the effect of if m = 1 ; this proves that the singularity of 

< 0, > is no worse than ^-h'+i-Sm,,!) _ Note that if T mj - is not singular (v = 0) then 
the same reasoning shows that < 0, h m j > will also be nonsingular. In summary, we have 
proven 


Proposition 6.2 If the asymptotic behavior ofT m j satisfies 


for an integer v > 0 , then 


with 


If v = 0 in ( |118|) , then 


lim Y\ T m J < oo 

7—>0+ 


lim y 1 " I < ib, h m ,■ > I < 00 
7^0+ 


v = 


v m = 1 
v + 1 m > 1 . 


lim | < - 0 , h m j > | < 00 . 

7—>0+ 


(118) 


(119) 


( 120 ) 


( 121 ) 


This result is useful because u' in 


often turns out to be less than the index of the 


integrand of < 0, h m j > and thus the Proposition can provide sharper estimates of the 
actual divergences. 


It is very important to recognize that the conclusions in (|119|) and (|121|) do not in 
general apply to < 0, h* m - >. Even though h m j and h* n i have the same index, the most 
singular terms in < 0, h m j > are not necessarily the most singular terms in < 0, /i* (J >. 
For < 0, h* m • >, if we wish to avoid explicitly evaluating the integral, our only means of 
estimating the divergence is by calculating the index of the integrand. 
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6.2 Calculation of the index 


These definitions allow the index of I m j and h m j to be calculated recursively to all orders. 
We illustrate this index calculation for the first few levels of Table I; our most general results 
are stated in Theorem [i . 1 1 and Proposition 0 - 

At the top of Table 1 we have the eigenfunctions defined by linear theory, these have a 
single pole in the upper half plane so 

Ind [ip{w)\ = Ind [V’(cu)] = 0 . ( 122 ) 

Next we consider J 2) o and /r 2j0 ; from (|79|) and (^3|) we calculate 

Ind [I 2 , oM] = 0 Ind [fi 2 ,oM] = 1- (123) 

The index of J 20 follows from (|122| ), while the expression for fi 2 0 has two terms; the first has 
index 1 and the second has index 0 , hence the index of h 2y o is 1 . 

Since fi 2j0 only has poles in the upper half plane, the integrals T 2 0 and < if), fi 2j0 > have 
no pinching singularities and remain finite as 7 —> 0 + : 

T 2 ,o ~ 7° <xp,h 2l 0 >~ 7°- (124) 

Note that the integrand in (|124|) has index equal to 2, but there is no divergence, this 
illustrates the second part of Proposition |6.2| . By contrast the integrand of < ip^ip* > has a 
smaller index than < ip, h 2 0 > but yields a singular integral ([85]) as 7 —> 0 + : 

< -0, ip* >~ 7 -1 . (125) 

The fifth order coefficient p 2 provides an instructive exercise in applying the index to 
estimate the singularity of a more complicated expression. The general form of the coefficient 
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(not assuming f 2 i = 0) is 


P 2 


-(27T iKl) [// < 0, /l 2) l > +T* j0 < 0, h 2 , 0 > + f r 2l r^o < 0, 4,0 > 

+f 2 l [ r 2 ,l < 0,0* > +r 2 ,0 < 0, 4,0 >] + /s« r 3,0 < 0, 4,0 >} 


(126) 


where 


4o = 27TiX(3/) [/; /i 2) o + / 2 *,0r 2)O ] (127) 

/ lj0 = (2mKl) P ± {/z 0 2 , o + f* 2l r r 2 , 0 } (128) 

0,i = 2pi/i 2 , 0 + 27riT7(2/) [/z/i 3 , 0 + /;*(/zi, 0 + Ti, o 0) + 03^3,00*] . (129) 

We can immediately estimate the integral < 0,4o > from the fact that the integrand 
has index 2 

< 44 , 0 >~ 7 - 2 . (130) 

The two terms in J 3 , 0 have indices of 1 and 0, respectively, and the poles of I 30 are all in 
the upper half plane. Writing out the expression in (|128|) 

Pj_ {04,o + / 2 /0* 10,0} = fih 2,0 + f 2 i 0* ^ 2,0 — 0 < 0, [/«4,o + f 2 i 0* 10,o] >, (131) 

we find terms of index 1, 0, and 1, respectively. The third term is already singular due to 
< 0,0* > and the second term has a pole in the lower half plane. These observations imply 


Ind [/ 3 ,oM] = 1 Ind [ZyoM] < 1, (132) 


and from (|75|) 


i,o ~ 7 


-2 


3,0 ~ 7 


0 


(133) 


The inequality in (|132|) reflects the fact that there are two terms of index 1 and we have 
not ruled out the possibility of a cancellation. The indices of h 30 and 4,o are now easily 
evaluated from ([76]) 


Ind [/i 3 , 0 (u;)] = 2 Ind [hi, 0 (u;)] < 2, (134) 
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and the corresponding integrals are estimated following (|119|) 


< 0, hz ,o > ~ 7° 


< 0 , hifi >~ 7 2 . 


(135) 


For p 2 we also need to estimate < 0 , h\ 0 > and this is done by noting from (|134|) that the 
integrand has an index of (at most) 3 so 

< h i,o >~ 7 " 3 - (136) 

The five terms in J 2 .i all have index 2 except T 3O 0* which has index 0, hence 

Ind < 2. (137) 

The terms involving h lj0 and 0* have poles in the lower half plane and therefore produce 
pinching singularities in r 2 p. By comparing (|75|) for T 2) i to (|125|) and (|135|) we can estimate 
a maximum singularity of 

r 2 ,i ~ 7 - 2 (138) 

which in turn implies Ind \h 2 ^{ui)] < 3, and 


< 0, h 2 ,i >~ 7 


-3 


(139) 


from (|76|) and ( |119|) , respectively. 

The result in (|139|) gives the divergence of the hrst term in p 2 , the asymptotic behavior 
of the remaining five terms can be evaluated from the foregoing estimates: 


< 0, ft-2,0 > 

~ y 2 

(140) 

< 0, h 3i o > 

^ 7 ° 

(141) 

L < 0,0* > 

~ 7 —3 

(142) 

< V’, hi ,0 > 

~ J~ 3 

(143) 

A 

to * 

o 

V 

~ 7 -2 . 

(144) 
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Hence barring some accidental cancellation among the most divergent terms we expect P 2 ~ 
7 -3 for the general case when /a 7 ^ 0 ; more precisely we have proved 


lim 7 3 |p 2 | < 00. (145) 

7^0+ 

6.3 The main result 

Our main result on the singularity structure of the amplitude expansions can now be proved. 
In the preceeding sections, we have calculated the indices of J mj - and h m j through the first 
two levels of Table I, and determined the singular behavior of the cubic and fifth order 
coefficients. This information is now inferred for the entire theory to all orders. 


Theorem 6.1 For j > 1, the singularities of the coefficients in the amplitude expansion 
satisfy 

lim 7 2 ! -1 \pj\ < 00. (146) 

7^0+ 

For m > 1 and j > 0, the indices of / mj and h m j obey 


Ind 

— Jm,j T 2<5 TO) 1. 

(147) 

Ind [hmj] 

^ Jm,j + 1 + 

(148) 


where J m j = m + 2 j — 2, and the integrals in © satisfy 


lim 

7^0+ 


lim 7 Jm ’ 

'+3<5 mj i i-p 

\ L m,j 

7 —> 0 + 


l - ! - 2^777,7 

' 0 , h m j 


< 00 

< 00. 


(149) 

(150) 


We are only able to prove the results in ( |147| ) - (|148|) as upper bounds on the index, but 
we expect them to hold as equalities in most cases. Note that the bounds in ( |149|) - (|150|) 
give weaker divergences than one would estimate simply from the index of the integrand as 
calculated using (|147|) - (|148| ). For example, from ( |148| ) the integrand of < -0, h m j > has 
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index J m j + 2 + 2<5 m> i, but the bound in ( |15C| ) assures us that the divergence of this integral 
will always be at least one power of 7 less. There is no contradiction since the index only 
gives a strict upper limit on divergence exhibited by a given integral. 

Proof. 


1. The relations (|146| ) 


have been checked explicitly for I 2 0 , h 20 , r 2j0 , 


and pi in the analysis of section |5.1| . In addition, in section |6.2| they have 
been vended for I lj0 , h lfi , T li0 , / 3 ,o, h 3j0 , T 3i0 , I 2} 1 , h 2 , 1 , T 2 ,i, and p 2 . This 
proves the theorem for the hrst two levels of Table I. 

2. We extend (|146|) - (|150|) to all pj, and { h m j , I m j} by induction using the re¬ 
cursion relations. Assume that ( |146| ) - (|150|) are true down to some arbitrary 
level of Table I, and consider what the recursion relations imply for the co¬ 
efficients evaluated at the next level. Let h m >^ and I m \j' denote coefficients 
at the next level that can be evaluated from lower order quantities which 
satisfy (|146|) - ( |150| ). We hrst consider the various recursion relations for 
h m >j> and and prove that if ( |146| ) - ( |15(J| ) hold on the right hand side of 
these relations, then ( |147|) - (|150|) will also hold for the coefficients obtained 
on the left. We organize this part of the proof by noting the importance of 
establishing two relations for I m <.p'. 

(a) the index identity (p~47|) 


Ind ^ Jm',j' T 2 S m ',!• 


(151) 


(b) the estimate 


hm r y'^m r , 1 

7 —> 0 + 




du- . 

-00 (CJ — Q,) 


< OO 


(152) 
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where a is any pole of the form described in ( | 107 | ) (in the upper half 
plane). Note that ( |151|) implies that the integrand in (|152| ) has an index 
°f + 2 8m ',1 + 1 ; thus ( | 152 | ) gives a divergence which is less by one 

power of (l/y) than would be “naively” predicted using the index of 
Irn' .j' in ( |lol| ) • 


If we can prove ( |151| ) and (|152|) for then the remaining properties 

( |148|) - (|150| ) are easily obtained as follows. First, (|149| ) follows from 


and the identity d75|). Next (|148|) follows by applying (|151| ) and (|149[) to the 
identity (|76]). Finally (|150| ) follows from Proposition [n2] and (|149|) . Thus 
the crux of the matter is to verify (|151| ) and ( |152|) ; this is done using the 


recursion relations for i.e. using (j77j) - (|8l|). 


3. Verification of ( p~5T|) . 

The index of each term on the right hand side of the recursion relations for 
I m r jr(uj) can be evaluated by applying (|146|) - (|150|) ; this exercise shows that 
all of these terms have an index less than or equal to + 2This 

establishes ( |151|) for I m /y. A few examples from o illustrate these index 
calculations. 


(a) The terms in i) j that depend on pj have the form [(2 + k)pj_ k + (1 + 
k)p*_ k ] hiy, from (|148| ) we have Ind [hi )k \ <2k + 2 and from ( |146| ) the 
singularity of pj-k is determined, hence 


Ind [ Pj _ k h ltk ] < 2{j -k)- 1 + (2k + 2) = 2j + 1 (153) 


which is consistent with ( |151|) for (■ m',j') = (1 ,j). 
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(154) 


(b) Next consider the terms 


j+2 j—m+2 

fml P± [^m-l,n^rn J-m-n+ 2 ] 

m=3 n =0 


which expand to 


— r m j_ Tn _ n _)_2 [h m _ ln < l,n V 7 ]- 


The index of each term can be evaluated 


Ind 


^-l/mj-m-n+2 — 2j 


Ind 


^ V’j hm—l,n > ^-m,j—m—n+2'P A 2j + 1; 


(155) 


(156) 

(157) 


thus Ind 


^ i\h* rn _i n T m ,j-m-n+ 2 \ < 2j + 1 which is consistent with 
( |151|) for ( m! ,j') = (1 ,j). In this example, both indicies are obtained 
from ( |148| ) and (|149| ). The divergence of T mj -_ m _ n+2 ~ f y-(2j-m-2n+2) 


h 


m—l,n 


< m + 2n — 2. 


is estimated from ( |149| ), and (|148|) gives Ind 
Hence the integrand of < ip, h* m _ l n > has an index no greater than 
m + 2n — 1, and we conclude that the worst possible divergence is 
< 7 _ ( m+2n_1 ). These estimates immediately imply ( |156| ) 

- (EH- 

(c) A more subtle example is the term /i 2 j in (|77|); from (|6T]) we have 
that the index of Pj_ h = h, 2 ,j — ip < ip, h, 2 j > is equal to the largest 
index obtained from h 2 j and < ip,h 2j - > '0( a O- These terms have 
indices 


Ind [h 2 j] < 2 j + 1 


Ind 


< ip, h 2 j > ip(u) < 2 j + 1, (158) 


so Ind [P_l h 2i j] < 2j + 1 which is consistent with ( pL51|) for ( m',j') = 
(1 ,j). The first index in (|158|) comes from (|148|) and the second follows 
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from ( |150| ) and ( |123| ). This example is delicate because ( |158|) requires 
that we use ( |150| ) to estimate < if, h 2 .j >~ q _ ( 2 -?' +1 ) rather than ( |148| ) 
as in the previous calculation ( |157| ). 

Proceeding in this way the index of every term in the recursion relations 
(0)-© can be calculated and (|151|) verified. 


4. Verification of 52|). 

We prove ( |152|) by integrating the recursion relations for and ex¬ 

amining the resulting recursion relations for 

(/((n) I m i jfiu) 


du- 


(159) 


(u — a) 

The right hand side of the recursion relations m-m involve products of 
if, if*, h m j, or h* rrh] with integrals whose singularities can be estimated by 


applying (|146|) - (|15Q) as illustrated above. Now when we divide by (a; — a) 
and integrate as in (|159|) we get additional integrals of the following types 


g{u)if(u) 


/ du V \ 

'—oo yOJ — CVJ 


g(u) 


du- 

I—co yUJ — (X.J 


du 


g(w) KAu) 


du 


g(u) 


(160) 


(161) 


(cn — ck) J —oo (cn — of) 

The first integral in ( 160 ) has no pinching singularity and the second has 
index 1 and diverges as 7 -1 ; the singularity of the integrals in ( 161 ) can be 
estimated from ( |148|) and (|150|) , respectively: 


lim W m J+ 2+2 ‘ 5 "M 

7—>0+ 

lim ryJmJ+l+ZSm^ 
7—>0+ 


roo 

' — oo 
roo 


du 


g(v) 


(u — a) 
Cj{u) hmji^u) 


du- 

l-oo (cn — a) 


< 00 


< 00 . 


(162) 

(163) 
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With these estimates the terms contributing to (|1 59| ) can be shown to di¬ 
verge at worst like 

-j • < 164 ) 

This establishes ( 11521) for and completes the induction argument that 

( |147p - (|150|) are propagated by the recursion relations. 

5. Now we finish the proof of the theorem by verifying that the estimates for 
the integrals in (|149|) and (|150|) imply the estimate of pj in (|146|) ; this is a 
straightforward application of (|149| ) and (|150|) to the right hand side of the 
recursion relations (TOi). A few examples from (|70[) suffice to illustrate how 
this is done. 


(a) Consider the term < >; from (|150| ) the singularity of this 

integral is at most 

J 2 , j —1+1 / 1 \ 2j — 1 


< ijj, h 2 ,j -1 > 


4 


-4 


(165) 


which is consistent with (|146|). This example is another case where the 


divergence of the integral is best characterized with ( 150 ) rather than 
inferred from the index of the integrand which is 2 j. 

(b) Other terms involve products of integrals, for example, T 2 ,j-\ < ip, ip* >■ 
Using (|149| ) and ( |125|) we determine that the singularity of this product 
is at most 


r 2j --i <ip,ip* > 


4 


J2.1— l+l 


(166) 


which is the same as in (|165|) . 

(c) The terms 


3 + 1 3 

m =3 

—m+1 

V r 

/ v 1 m,n 
72=0 
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are a final example. Applying 
gularity of 

r m ,n < V’) lj-n-m+l ^ 

which is also consistent with 


and (|149|) shows a maximum sin- 


^ \ ^777, TlS - J TTl — 1, J — 71 — 771-(- 1 ~t - 2 / | \ 2 j 1 


.7, 


.7, 


(168) 


In this manner every term in the recursion relation for pj is shown to satisfy 


the general estimate in ( |146|) . 


□ 


6.4 Implications of Theorem fSTTj 

Our discussion in section |5.1| of the significance of the coefficient singularities did not consider 
terms in the amplitude equation beyond third order. If all such terms are included then the 
amplitude equation in polar variables (|90|) becomes 


A = 7 A + 

3 = 1 
oo 

-0 = Q — ^2 Im (pj)A 2 - 7 . 

3 = 1 


(169) 

(170) 


With the asymptotic behavior of the coefficients is given by Theorem |0 

[ a j + ^(7)] 


Pj = 


7 


2j—1 


(171) 


we re-express the dynamics in terms of the rescaled amplitude (@) A(t) = y /3 r( 7 ' 5 t) 

^ = 7 1 - 5 r + £ 7 2 ^- 1 )- 5+1 Re(a,)[l + 0(7)]r 2 ^ 1 
dr j=i 


di[> 

dt 


= — y~) 7 2j ^ ^ +1 Im(aj)[l + (9( 7 )]r 2j . 


(172) 

(173) 


1=1 
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Now we see that the singularities of pj are such that the choices <5 = 1 and (3 — 1 remove 
them to all orders, leaving an asymptotically finite system 


Hr 00 

— = r + Re(a ,-)[1 + 0{^)\r 2j+l 

dT j=i 

dip 

dt 


= — 7 ^2 + G>(l)}r 2 j■ 


(174) 

(175) 


3 =1 


In ( 174 ) there is a further remarkable consequence of the singularities, namely that the 
terms at higher order in r no longer appear to be negligible even for small 7 . This conclusion 
is provisional in the sense that we have not proven aj 7 ^ 0 at high order. Nevertheless, our 
analysis of p\ and P 2 shows that in general cq and <22 are non-zero, and inspection of the 
recursion relations shows no evidence that the series will truncate. This feature does not 
affect the phase evolution (|175|) quite so dramatically. If we assume a time-independent 
solution r = R rw of ( 174 ) then the nonlinear frequency of this state 

OO 

n™ = n - 7 £ M«;)[1 + o(i)\R 2 r i (176) 

3 = 1 

is very close to the linear frequency for small 7 . Such solutions are called “rotating wave” 


states as they describe a solution 




7 


Rr 




T + a;) + ^ Rl w e - 1 + ( u, y 2/3 R 2 rw ) e 


iW 


(177) 


£7 


h m (u,y w R 2 r J e‘ 


m =2 


+ c.c. 


that appears time independent in a rotating frame O' = 6 — Q rw t/l. 
The asymptotic behavior of T mj - is given by Theorem |6.1| as 

[bm,j + ^( 7 )] 


r . — 

L m,j j . 1 ox . 


(178) 


where the constants b m j are independent of 7 . The expansion of T m thus has the asymptotic 
form 

OO / 1 \ m 2 + 3577^,1 OO 

(179) 




£[ 6 mj + 0{n)]r 2 K 
3=0 
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This motivates the definition of the nonsingular coefficient T m 


m-2+3Sm t i 


r m fcr) = 


. 7 , 


f m (r 2 ); 


(180) 


T m will play a role shortly in our discussion of Daido’s order parameter. 


7 Special Cases: fa / = 0 


The general analysis shows that if j' 21 , 7 ^ 0, then the resulting singularities require a rescaling 
of the amplitude by a factor of 7 . In this section, we consider some of the possiblities that 
can occur for couplings with j '21 = 0. First, we note that the general recursion relations 


simplify somewhat when f 2 i = 0 ; the terms proportional to j '21 are explicitly shown in ([70]) 
and ( |77| ) - (| 8 TD, and we simply drop them. 


7.1 Absence of the cubic singularity: / 2 Z = 0 and fa ^ 0 

When f 2 i = 0 so that the cubic coefficient is non-singular, then if f 3 i 7 ^ 0 the first singularity 
appears in the fifth order coefficient p 2 ~ y~ 2 . The general structure of the amplitude 
expansion for this situation is characterized by the following proposition. 


Proposition 7.1 Suppose that f 2 i = 0 and f 3 i 7 ^ 0, then for j > 1, the singularities of the 
coefficients in the amplitude expansion satisfy 


lim y 2 - 7 ' 2 \pj\ < 00 . 
7^0+ 

For m > 1 and j > 0, the indices of I m j and h m j obey 


Ind [Ai.j] 77 Jm,j T 2<5 mil 
Ind ^ Jrn,j T 1 T 


(181) 


(182) 

(183) 
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where J m} j = m + 2 j — 2, and the integrals in © satisfy 


lim 5 m, 26 j , 0 ) | F | < qq 

7—>0+ 


lim Y m ’ j+2Srn ’ 1 < ' 0 , /i m j > 

7—>0+ 


< OO. 


(184) 

(185) 


Compared to the conclusions of Theorem | 6 .1| the bounds on the indices are the same so there 
is nothing to prove. The estimates of the singularities in ( |181|) and (|184j ) - ( |185[) are softened 


by one factor of 7 . Since J 2j o = 0, there is no singularity in T 2j o even in the general case 
and the odd correction term <5 mi 2 <5y 0 hi (|184j) adjusts the exponent to allow for this. In fact, 
r m , 0 , m > 1 is nonsingular in general, although we have not used this. 

Proof. 


1. The results in (|181|) - (|185|) have been established by the calculations in 
sections IOI and fv 2 |: they are extended to all orders by induction. The 


argument is entirely analogous to the proof given for Theorem ^T| and we 
indicate only the necessary modifications and omit the details. 

2. The only adjustment required is that the second crucial identity (|152|) in 
the general proof is now replaced by the estimate 


lim ,j' + 2 S m' ,l- 0 -- S m' ,2 S j' ,0) 

7—»o+ 


r °° 

—00 (to — ck) 


< (X) 


(186) 


where a is any pole of the form described in (|107|) (in the upper half plane). 
The first index identity (|151|) is still appropriate: 


Ind Jm\j' T 2 


(187) 


With this adjustment the repetition of the previous argument is straight¬ 
forward. 


□ 
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7.2 Implications for scaling 


The consequences of this Theorem for the scaling of the amplitude are ambiguous. If we 
rewrite the radial equation from (|172|) with 5 = 1 

dr 


dr 


= r 


1 Re(p J )r 


2j+l 


(188) 


3 = 1 


and then insert the estimated singularity pj ~ cq-q 2 ^ +2 we find 

dr 


dr 


= r + 1 - )+1 Re(aj)r 2 - 7+1 . 


(189) 


3 =1 


Now the requirement that all the nonlinear terms should be nonsingular, 2 j{(3 — 1) + 1 > 0, 
results in a j-dependent estimate on the exponent: (3 > 1 — 1/2 j. We know the fifth order 
singularity occurs so from j = 2 we have (3 = 3/4 as discussed previously in Section fjT2| . 
However our analysis of the recursion relations is not sharp enough to determine if the leading 
singularities are actually present at every order. Consequently, we have a range of possible 
values of (3 


<(3< 1. 


(190) 


This range can be reduced to (3 = 1 if we require that the rescaling also yield a finite 
expansion for T with coefficient singularities as estimated in (|184| ) . However it must be borne 
in mind that these singularities are upper bounds. Our firmest prediction remains of the 
range of values in (|190| ) since the bound (3 > 3/4 is conclusively required by our fifth order 
calculations and the the bound (3 < 1 follows from our analysis of the general case. 


7.3 Single component couplings 

An additional case of interest is the circumstance which holds for the Kuramoto model 
/((/>) = sin</>; a coupling that lacks all harmonics of the Fourier component fi driving the 
instability. In this exceptional case the amplitude expansion is entirely free of singularities. 


47 








Proposition 7.2 If f n i = 0 forn > 1, then the coefficients pj in © are all nonsingular: 

lim \pj\ < oo (191) 

7 ^ 0 + 

for j > 0. 

Proof. 


The functions I m j(u) and h mj (uj) have poles in the upper half plane only and 
consequently the integrals appearing in the coefficients pj are all free of pinching 
singularities in the limit 7 —> 0 + . Also T mj - is nonsingular from (|75|) in light of 
the discussion following (|96|) . Hence there are no singularities. 


□ 


The absence of singularities results immediately in the conclusion (3 = 1/2 for the scaling 
exponent. 


8 Evaluation of Daido’s order function 


Daido’s definition of the order function H(9) is a direct generalization of the familiar order 
parameter of Kuramoto and is motivated by the mean field form of the phase dynamics 
®.0 |15| , fj4|| He assumes there is a component of the population collectively entrained 
with frequency Vl e and introduces the shifted phases fj = 6j — Q e t, then phase equations (H) 
can be rewritten 

ifi = Ai - H(i(ji,t) +£i(t) (192) 


where Aj = u;* — Fl e and H(9, t) is the mean field 


H(9,t) = - f” e 

n=—o o 


—inO 


N 




N 


3= 1 


(193) 
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The additive noise term £j(£) in ([l)is omitted in the analysis, and for large N, Daido assumes 
Vt e can be chosen so that the limit 


H(9) = lim H(9,t) 

t—>00 


(194) 


exists. The time-asymptotic mean field then defines the order function H{9). A norm of H , 

r2n 

\\H\\ 2 = / d9H{9) 2 /2ir, (195) 

Jo 

serves as an order parameter for models with multi-component couplings. 

It is interesting to evaluate H using our continuum description. The average over the 


population in (|193|) can be expressed in terms of p(9,LU,t)\ 


1 N /*27T f'OO 

lim y e » (s i (() -M= d9 dcu^(a;)p(0,u;,f)e in(0 - nef) , 
N-^oo N “[ Jo J- oo 


(196) 


and this yields a general expression for the order function 


°° p2tt roc 

H(9) = — lim y f n e~ ine d9' du'g{u')p{9', u/, t)e in ^ 6 . (197) 

J ° J 


Let H u denote the order function obtained from solutions on the unstable manifold (|6a|); a 
straightforward substitution from (^) then yields 


H u (9)/2it = -Urn {fie 




a(t)* + a(t)*a(t) / dcu g{uj)h\{uj , cr(£))’ 1 


(198) 


y frme 


-iml( 9 +n e t) \ m 


du g(u)h m (u, a(t))*) + cc > ; 


m=2 


note that we have set fo = 0 for simplicity. 


Daido’s assumption of time-independent behavior (|194|) is naturally satisfied when the 
system is described by a “rotating wave” solution of the form (cf. 


otrw{t) = R rw e l<t>0 e 


% ruit 


(199) 
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with R rw and Q rw the time-independent amplitude and frequency, respectively. For these 
asymptotic states by choosing Vt e = Q rw /l, the limit in (|198|) exists and the order function is 


m 


2t r 


ft 1 l>R rw e- i *°(l + 1 2l> R 2 rw r 1 )e ,w +'£ f' m , e"'”*" T m e‘ 


iml6 


m =2 


+ cc (200) 


with norm given by 


| Hf = (2tt) 2 7 2/3 # 


>2 

l rw 


fil 2 (i + 7 2 "«L r i) 2 + E l/ m il 2 7 (2 ” _2) ' 5 ^, 2 ” _2 ri 


( 201 ) 


m=2 

Our results on the possible values of (3 imply corresponding predictions for the scaling of 
11//’ll near onset. There are several different cases depending on the coupling and whether or 
not noise is included; a summary appears in Table II. 


1. Suppose that either D > 0, or D = 0 and we have a “single component” coupling: 
fi 7^ 0 and f m i = 0 for m > 1. Then as 7 —> 0 + , there are no singularities and (|201| ) 
gives 


\H\\ = 27T7 PR r 


1 + 0( 7 2/3 )l ~ (K - K c f 


( 202 ) 


with (3 — 1/2. In each of these circumstances, the prediction \\H\\ (/I - A c 7 2 is in 

agreement with Dai do’s analysis of ||ZZ||.||1 5|| 


2. Suppose that D = 0 and we have the general case with fi ^ 0 and / 2 / 7^ 0. The 
coefficients of the amplitude expansion have the generic singularities (|146|) and the 
coefficients T m are also singular (|180|). Substituting (|180|) into (|201|) we hnd 


\Hf = (27r) 2 7 2/3 i? 


2 

rw 


f i i 2 (i+7 2( ' , - 1, <fi«)) 2 + E i/ m ii 2 7 2,m - 1)(,5 - 1, ^r 2 f; 2 .. 


m =2 


in this case we have determined that (3 = 1 which yields 


| H\\ = 27T7 R r 


f/| 2 (l + OM 2 +E \fml\ 2 R 2 r™~ 2 T 2 


m =2 


1/2 


(203) 


{K-IQ 1 . (204) 
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This is the scaling for \\H\\ Daido found for transitions at l = 1 when there was no 
noise and the coupling had components at l = 2. We obtain this result for a transition 
at arbitrary l if the coupling satisfies / 2 ; 7 ^ 0 and D = 0. 


3. Suppose that D = 0 and we have the special case with fi 7 ^ 0, / 2 z = 0, and 7 ^ 0. In 
this case, our analysis of the recursion relations leaves considerable uncertainty in the 
value of f3 as indicated by the range (|190|) . 


4. Suppose that D = 0 and we have the special case with 0(2) symmetry and also a 
coupling that satisfies /(</>) = —/(</> + 7 r/Z). This requires f n = 0 unless n is an odd 
multiple of l and the non-zero components must be imaginary. I 11 this case again one 
has typically // 7 ^ 0 and fsi 7 ^ 0 , and our analysis of the recursion relations specifies 
/3 only within the range (|190|). However this range nevertheless excludes the value 


/3=1/2 Daido obtains for such symmetric systems. |p~3|. [ 34 ]. | 35 | The reason for this 
apparent discrepancy has not been found. 


9 Discussion 


Our analysis shares a key conclusion with Daido’s study of scaling based on the order func¬ 
tion. I 11 the absence of noise, the presence of a coupling component f 2 i at the first harmonic 
of the critical mode number l will slow the onset of synchronized behavior so that the phase- 
locked component of the population scales like (K — K c ) rather than (K — K c ) 1 / 2 as in the 
single-component coupling models. However, both approaches suffer a common deficiency; 
neither provides a simple explanation of why the onset of synchronization is slower if the 
coupling has a second harmonic component (/ 2 Z 7 ^ 0). This shift seems to be a true many- 
body effect requiring large N, but a dynamical mechanism has not been identified. Perhaps 
further numerical studies will shed some light on this issue. 


51 











For a general coupling with f '21 7 ^ 0, our calculations show two effects of adding noise to 
the phase dynamics. The diffusion term tends to suppress the linear instability so that K c 
increases, but in addition the noise prevents the poles from actually reaching the real axis 
and therefore removes the coefficient singularities; this modifies the scaling (K — K c ) ,a of the 
synchronized state changing (3 from 1 to 1/2. The interaction of these two effects raises the 
intriguing possiblity that under some circumstances, adding noise could actually yield an 
enhancement in the level of synchronization. In some recent studies of mean field oscillator 


models it has been noted that noise can enchance a collective oscillation. 36 


The linear stability analysis of section |2] clearly shows that for multi-component couplings 
one can expect multiple linear instabilities and this will allow for studies of the interaction 
between two different synchronizing transitions. In other settings the study of simultaneous 
linear instabilities, or codimension-two mode interactions, has been one of the keys to ob¬ 
taining an analytical theory capable of describing the system beyond the threshold of the 
first instability. Such a theory can identify possibilities of secondary bifurcation and predict 
the appearance of states with nontrivial time-dependence. Despite the widespread interest in 
synchronization, and the large literature on the Kuramoto model, there has been no analysis 
of the codimension-two bifurcations that arise when more than one mode of synchronization 
is possible. 

Although the justification for a phase model varies from one context to another, in the 
setting of coupled oscillators it requires two essential assumptions: weak coupling and a 
frequency spread g(u) that is “sufficiently” small. |l], §, |3| Moving beyond these assump¬ 
tions necessitates restoring the dynamics of the oscillator amplitudes rj(f), and this adds an 
additional term to the kinetic equation of the schematic form 

d 


dr 


(rp) 


(205) 
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where r = ( 77 ,... ,rjy) are the oscillator amplitudes and r is the non-stochastic component 
of the amplitude dynamics. In this situation, p describes the joint distribution of oscillator 
phases and amplitudes. The effect of such a term on the synchronizing transition can be 
studied by analyzing the effect on the singularities of the nonlinear coefficients. 
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A Proof of Proposition 1671] 


Let zq = (f2 T i^)/l denote the root in ( |107|) - (|108|) . If run = 0 then there is no pinching 
singularity and the integral in ( |111|) has a finite limit; hence the limit is zero if J > 0. 

Assume mn > 0, then for integrals with m + n = 2 the limit can be evaluated from the 
Plernej formula Jl7|| after a partial fraction expansion. The expansion isolates the singularity 
at 7 = 0 , 

o ✓"/ / $ / /it n I t /* no i ^A //I i rvi si / i /-A if it 

(206) 


f°° du p(u) 

—i/l 

f°° du 4>(ui) 

f°° du (j>(u) 

J—oo (c<J — Q?i)(cd — fil) 

+ Ci + 2/1 

J—oo (id — Q^i ) 

J- 00 (u-p{)_ 


and the limit (| 1 1 1 | ) yields 
° du <f>(u) 


lim 7 

7—>o+ 


/—00 {u — Qti){u — Pi) 


(5 1 + Ci + 2/1) 7 - 0 + 
27r\(p(Q/l — Kf 0 ) | 


r°° du p(u) 

I —00 ((U — Oil ) 


du p{u) 


/-00 (u — PI ) 

(207) 


($i + Ci + 2/1) 

The partial fraction expansion in (|206|) expresses an integrand of index 1 in terms of 


(208) 


integrands of index 0. More generally, for integrands 

roc 

/ dup{u) Q{u) 
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with index m + n — 1 > 1 , by expanding the integrand in partial fractions, they can be 
re-expressed as a sum of two integrands with index m + n — 2 that are multiplied by an 
overall factor of 7 -1 . Thus, if the partial fraction integrands of index m + n — 2 satisfy ( | 111 | ) 
with J — m + n — 2, then the original integrand will satisfy ( |111|) with J — m + n — 1; the 
value of J must be incremented by one to allow for the factor of 7 -1 . 


A simple induction argument along these lines establishes ( 111 ) for all rn + n > 1. In the 


limit in (Ml) is found to be non-zero in general so, barring an accidental cancellation, if J 
is replaced by J — 1 then the modified limit will diverge as 7 -1 . 

B Derivation of the model 


Basically, the kinetic equation arises as the first member of a coupled hierarchy of equations 
quite similar to the BBGKY hierarchy well known in the kinetic theory of gases and this is 
probably the most elementary systematic approach. There are alternative derivations using 
path integral methods. pBf 

The solution of the system of N stochastic differential equations in (|lj) is a Markov process 
with transition probability V(9, to, t\9 u/, t') where 6 = (0 1; ..., 0 N ) and u = (uq,..., cyy). [j39j 
This transition probability satisfies a linear Fokker-Planck equation 


dv 9 (ev) = d 9 


dt 39 


39 39 


(209) 


subject to the initial condition 


V{9,uj 1 t'\9\uj' 1 t’) = 5 N (9 — 9') 5 N (uj — u'). 


( 210 ) 


The notation 9 = [9 1 , ..., On) in ( |209| ) refers only to the non-stochastic part of the dynamics 

K N 

e l = u Ji + -J 2 f(o J - 9 l ). ( 211 ) 

iV 3 =1 
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Given an initial normalized ensemble of populations So (9,uj), the ensemble for t > to is 
determined by 

£(0,u,t) = J dJ J d9' V(9,uj,t\9',uj',to) £ o (0',a/); (212) 

equivalently, this is can be written as an evolution equation for S(9, u>,t) 


dS_ 

~dt 



de de 


(213) 


with initial condition £ 0 . The original phase dynamics (|lj) is unchanged if we exchange a 
pair of oscillators ( 9i,Ui ) <->• ( 9j,Uj ) and it is natural to select initial ensembles So with this 
same invariance under pairwise exchange. In addition, for our problem, the values of the 
frequencies a;* are assumed to obey a fixed distribution g(u). This requirement necessitates 
a suitable choice of initial ensemble: 

N 

£o(0 > u ) — Go(9, u) g(ujj) (214) 

3 =1 

where the normalization of Sq(8,uj) implies 


1 = J de G 0 (9,u). 


(215) 


A reduced distribution p s for the ensemble is defined in the usual way 


p s (0i ,..., = j d0 3+ dO N j <Li s+1 --- du N £(0,u,t), 


(216) 


and one can derive a hierarchy of coupled equations expressing the evolution of p s in terms 
of p s and p s+ 1 by performing the integration in ( 216 ) on the evolution equation ( |213| ). For 
our purposes it is sufficient to find the equation for pp integrating ( |213| ) over (0 2 ,... ,0n) 
and (u3 2 , ■ ■ ■ ,ujn) yields 


dpi d 

~dt + W 1 


(uq + Kf(0)/N)p 1 + — ~~~r ^ ( dd 2 [ du 2 f{9 2 - 9i)p 2 {9i,0 2 ,uJi,U2,t) 


N 


= B W (217) 
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In terms of the two-oscillator correlation function, 


62, Uli,Ul2, t) = Pl(#l, Wi, t) Pl(02,U)2, t) — P2 ($1; 82, Wl, U>2, t ), 


this becomes 


dpi d 

~dt + W 1 


{ui + Kf( 0 )/N)pi -\ - 1 --—-— j (I 62 J dU 2 f (02 — 9i) Pi(02) t) 


D^ + — 
def dOi 


K{N - 1) 
N 


dd 2 / du2f(e 2 — 9i)C(9i, 6 * 2 , uq, uq, t) 


(218) 


We obtain an autonomous equation for pi by letting N —> 00 and discarding the corre¬ 
lation term 


dpi d 

~dt + W 1 


UJiPi + Kpi J d02 J dU2f{02 — 9i) Pi( 02 ) ^ 2 , t) 


= (219) 


In some cases, one can rigorously prove that with a mean held interaction the effects of 
correlations are negligible in the limit N —> oo.|M| [41], ^]) Whether this kind of result holds 
for the models considered here appears to be an unsettled issue. Finally, the ensemble choice 
in (|214|) implies that 


r 2n 


d9 1 p 1 (9 u u 1 ,t 0 ) = g(uji), 


( 220 ) 


so it is natural to define p(9i,ui,t) by 


Pi( 0 i,wi,t) = gMp( 9 i,u u t). 


( 221 ) 


Substituting this expression for p\ into ( |219| ) yields the kinetic equation (p) - (|4|). 

It is interesting to note that the choice of initial ensemble fl214|) can be adapted to treat 


other forms of quenched randomness. [^] For example, we can allow for random variation in 
the coupling strength or even in the sign of the coupling between different pairs of phases. 
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Table 1: Order of calculation of h m j(u) and pj from p(u). The flow of calculation of the 
h m j(uj) is indicated by moving downward; the calculation of the p 3 proceeds by moving down¬ 
ward and to the left. For example, from p(u), /r 2 ,o can be calculated and then pi determined; 
hi 0 and /r 3i0 are calculated next from h 2y o and then /i 2 ,i can be evaluated from {/?i 0 , ^ 3 , 0 ? Pi}- 
This then determines p 2 , and so forth. For N > 2, pjy requires prior calculation of h m j for 
1 < m < N + 1 and 0<j<N — m + l — 25 m p. 



rn — 1 

rn — 2 

rn — 3 

rn — 4 

rn — 5 rn — 6 • • • 

Po 

V’M 

Pi 

- 

h- 2 fi 




P 2 

^ 1,0 

h 2 .1 

^ 3,0 



P 3 

hi,i 

h 2 p 

^ 3,1 

^ 4,0 


Pi 


h 2 p 

^ 3,2 

^ 4,1 

^ 5,0 

Pb 

hip 

h 2 p 

^ 3,3 

hip 

^ 6,0 

hr,.\ 




Table 2: Comparison of the exponent P obtain from the amplitude equation with the expo¬ 
nent (3d found from Daido’s order function. We assume the synchronizing transition occurs 
at mode number l for a coupling with // p 0. The bifurcation has 0(2) symmetry when the 
coupling and frequency distribution satisfy /(</>) = —/(—</>) and g(u) = g(—u), respectively. 


Noise 

Symmetry 

Coupling 

P 

t-H 

'■—o 

Comment 

D > 0 

SO(2) or 0(2) 

arbitrary 

1/2 

. 

generic case with noise 

D = 0 

SO(2) or 0(2) 

./21 ^ 0 

1 

1 

generic case without noise 

D = 0 

SO(2) or 0(2) 

f'21, = 0 , f:u 7^ 0 

[3/4,1] 

. 


D = 0 

0(2) 

/(</>) = -/(</> + 7 T/ 1 ) 

[3/4,1] 

1/2 

P and Pr, differ 

D = 0 

SO(2) or 0(2) 

fml = 0, m > 1 

1/2 

1/2 

single component case 
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